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The electrostatic interaction in ionic fluids is well-known to give rise to a characteristic phase be¬ 
havior and structure. Sometimes its long range is proposed to single out the electrostatic potential 
over other interactions with shorter ranges. Here the importance of the range for the phase behavior 
and the structure of ionic fluids is investigated by means of grandcanonical Monte Carlo simulations 
of the lattice restricted primitive model (LRPM). The long-ranged electrostatic interaction is com¬ 
pared to various types of short-ranged potentials obtained by sharp and/or smooth cut-off schemes. 
Sharply cut off electrostatic potentials are found to lead to a strong dependence of the phase behavior 
and the structure on the cut-off radius. However, when combined with a suitable additional smooth 
cut-off, the short-ranged LRPM is found to exhibit quantitatively the same phase behavior and 
structure as the conventional long-ranged LRPM. Moreover, the Stillinger-Lovett perfect screening 
property, which is well-known to be generated by the long-ranged electrostatic potential, is also 
fulfilled by short-ranged LRPMs with smooth cut-offs. By showing that the characteristic phase 
behavior and structure of ionic fluids can also be found in systems with short-ranged potentials, one 
can conclude that the decisive property of the electrostatic potential in ionic fluids is not the long 
range but rather the valency dependence. 


I. INTRODUCTION 

Over the last two decades, the scientific interest in 
organic salts with melting points near room tempera¬ 
ture, so-called room temperature ionic liquids (RTILs), 
has been growing steadily. Features such as a negligible 
vapor pressure and a remarkable thermal stability [1- 
4] promise future applications as solvents for syntheses, 
electrolytes in fuel cells, solar cells, and batteries, as well 
as in biomass processing [5-10]. Moreover, due to the tiny 
vapor pressure many applications are conceivable under 
ultrahigh vacuum [1, 3, 11, 12]. However, the technolog¬ 
ical use of RTILs requires an in-depth understanding of 
ionic systems. 

The Coulomb interaction, underlying these systems, 
acts repulsively for equally-charged and attractively for 
oppositely-charged ions and it decays oc 1/r, i.e., it is 
long-ranged. On the other hand, the pair distribution 
functions of an ionic fluid decay exponentially, which 
is called the Stillinger-Lovett perfect screening property 
[13-16] and which is a necessary consequence of the long 
range of the Coulomb potential [17, 18[. It is a challenge 
for more than a century now to develop quantitatively 
reliable theoretical descriptions of this peculiar combina¬ 
tion of properties. For dilute electrolyte solutions Debye- 
Hiickel theory [I9-21[ is typically a good starting point, 
whereas for dense ionic systems ion-pairing is relevant 
[13, 14[. From experimental work [22-25] as well as con¬ 
tinuum simulations [26-29], there is evidence that in par¬ 
ticular the critical behavior of ionic systems is very sim¬ 
ilar to the Ising universality class, which typically ap¬ 
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plies to systems with short-ranged interactions. This re¬ 
sult suggests that the long-range character of the bare 
Coulomb potential might be of minor importance for the 
critical properties of ionic fluids. 

The aim of the present work is to investigate the rel¬ 
evance of the long range of the Coulomb interaction 
for the whole phase diagram and the bulk structure. 
This is done by studying the lattice restricted primitive 
model (LRPM) for the Coulomb potential being trun¬ 
cated smoothly on a length scale 1/a and sharply at a 
cut-off radius rcut- The technical details of the model and 
the approach are presented in Sec. H. The long-ranged 
LRPM is well-known to exhibit tricritical behavior at the 
crossover from a first- to a second-order phase transition 
between a charge-ordered and a charge-disordered phase 
[30, 31[. The (L)RPM renders the ions as homogeneous 
charged hard spheres, which is of course not valid for ionic 
liquids in general, since they may exhibit charge as well as 
shape anisotropy. However, this work is concerned with 
the influence of the long range and the valency depen¬ 
dence of the electrostatic interaction, being omnipresent 
in all kinds of ionic systems. In particular, the impact of 
the long range will not depend on rather short-range fea¬ 
tures like charge or shape anisotropy, therefore we con¬ 
sider this simple model of ionic fluids as sufficient for 
the present investigation. It is shown in Sec. HI that, 
upon varying the smooth cut-off decay constant a and 
the sharp cut-off radius rcut, the short-ranged LRPM can 
be tuned from a system void of a charged-ordered phase 
via one exhibiting charge-ordered and charge-disordered 
phases but differing quantitatively from the long-ranged 
LRPM to a model with short-ranged interactions whose 
phase behavior and structure is quantitatively the same 
as for the long-ranged Coulomb interaction. Moreover, 
even the Stillinger-Lovett perfect screening property can 
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be fulfilled for suitable smooth cut-off potentials, an ob¬ 
servation which is not trivial in the context of short- 
ranged interactions. By showing that the characteristic 
phase behavior and structure of ionic fluids can also be 
found in systems with short-ranged potentials, it is con¬ 
cluded in Sec. IV that the decisive property of the elec¬ 
trostatic potential in ionic fluids is not the long range but 
rather the valency dependence. 


are equivalent in finite dimensions, the power law oc l/r 
is asymptotically preserved irrespective of the choice of 
the norm. The parameters p and q in Eq. (1) describe 
the norms to be used for measuring the distance deter¬ 
mining the interaction potential and the sharp cut-off, 
respectively. 


B. Grandcanonical Monte Carlo simulations 


II. MODEL AND METHOD 


A. LRPM with truncated Coulomb potential 


Consider the lattice restricted primitive model 
(LRPM) of univalent cations (valency Zq = -fl) and an¬ 
ions (valency Zq = — 1 ) with hard cores of diameter a 
occupying but not necessarily exhausting the sites of a 
three-dimensional simple cubic lattice with lattice con¬ 
stant a. Global charge neutrality is guaranteed by an 
equal number of cations and anions in the system, and the 
hard cores ensure that each site is either empty or singly 
occupied. Within the original LRPM two ions of species 
i,j € {©,©} at a distance r interact, besides the hard 
core exclusion, via the infinitely-ranged Coulomb poten¬ 


tial 
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where lb = 


is the Bjerrum length with 
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the electronic permittivity e, the Boltzmann constant k, 

and the temperature T, and where r is measured with 
the Euclidean metric. 

In the present work the relevance of the long-range 
character of the Coulomb interaction for the properties 
of ionic systems is studied by considering the implications 
of replacing the infinitely-ranged Coulomb potential by 
the truncated Coulomb-like potential 
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with the inverse temperature /? = I/(/cT), the decay con¬ 
stant a > 0, and the cut-off radius rcut- By means of 
Eq. (1) two different methods of cutting off the Coulomb 
potential can be studied: On the one hand, the fac¬ 
tor erfc(a||r||p), leads to a smooth cutting off on the 
length scale 1 /a, and, on the other hand, a sharp cut¬ 
ting off at radius ||r||q = rcut can be considered. In¬ 
spired by the short-ranged potentials appearing within 
the Ewald method [32], here the smoothening factors are 
chosen to be complementary error functions erfc(x) = 
2 / 7 r drexp(—r^). Obviously, the infinitely-ranged 
Coulomb potential corresponds to a = 0,rcut = oo. Be¬ 
sides the Euclidean norm ||r ||2 = + y'^ + z'^ to mea¬ 

sure distances, the 1 -norm ||r||i = |a:| + |j/| + \z\ and 
the supremum norm ||r||oo = max(|x|, jj/j, jzl) are consid¬ 
ered, which are more adapted to a lattice model since, 
for lattice vectors r, they lead to values which are inte¬ 
ger multiples of the lattice constant cr. Since all norms 


In order to discuss the thermal and structural proper¬ 
ties of the LRPM with the truncated Coulomb-like in¬ 
teraction Eq. (1), the packing fraction rj and the pair 
distribution functions gij{r), i,j G {©,©}, are deter¬ 
mined for cubic boxes V := {0, cr,..., (L — l)cr}^, L G N, 
with periodic boundary conditions using grandcanonical 
Monte Carlo simulations. The set of all configurations 
C of cations and anions occupying V can be expressed 
as the set of all maps ( : V —)■ {0, Z 0 ,Zq}, which re¬ 
sult in the charge /'(r) G {O,Z 0 ,Z 0 } located at position 
r G V, i.e., C(r) = 0 iff site r is empty, (((f) = ^0 iff 
site r is occupied by a cation, and /'(r) = Zq iff site r is 
occupied by an anion. Standard Metropolis importance 
sampling of the grandcanonical Boltzmann distribution 
P{Q oc exp(/3/iiV[/'] — 13H[Q) with the chemical potential 
/i, the total number of ions N[C,] = iV 0 [(/] + iV 0 [(/], and 
the Hamiltonian 

I3H[C] = ^ ( 2 ) 

r.r'eV 

r:^r^ 


on the set of all charge-neutral configurations C is ap¬ 
plied. Charge neutrality is preserved during the simula¬ 
tion runs due to insertions and removals of only neutral 
pairs of cations and anions, i.e., the number of cations 
Nq [C] equals the number of anions Nq [C] . 

The packing fraction is the average 


and the pair distribution functions are given by 




[c])) 
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with i,j G {©,©}, r being a distance measured in the 
1 -norm, I 4 heii(?') = 4(r/(T)^ + 2 representing the num¬ 
ber of all sites in the shell of 1 -norm distance r around 
a site, and IVij (r, [(/]) denoting the total number of all 
ordered pairs (ri,rj) of positions ri,rj G V being sepa¬ 
rated by a 1 -norm distance r = ||ri — rj||i and such that 
an ion of species i is located at position and an ion 
of species j is located at position r^ . Here the pair dis¬ 
tribution functions gij (r) are defined in terms of 1 -norm 
distances r because, for the present simple-cubic lattice 
geometry, this choice is most convenient in order to dis¬ 
tinguish charge-ordered and charge-disordered phases. 
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FIG. 1. Isotherms T* G {1.25,1.43,1.67} for L = 12, a = 
0,rcut/a = 1 using the 1-norm {jp = q = 1) m Eq. (2). For 
T* = 1.25 a distinct discontinuity of the packing fraction q 
at chemical potential p* « —2.3 can be observed, while for 
T* = 1.67 a continuous curve occurs. Therefore, the critical 
temperature T*, where the first-order gas/liquid phase tran¬ 
sition associated with the discontinuity of the packing frac¬ 
tion q terminates, is located within this temperature range: 
T* € [1.25, 1.67]. For each isotherm two sets of data points 
are shown, which correspond to either an initially empty sys¬ 
tem (large dots connected by solid lines) or an initially totally 
filled system (small dots connected by dashed lines). The in¬ 
set sketches the procedure to estimate the packing fractions 
at coexistence (see Sec. IIC). 


The equilibrium state of the LRPM is determined by 
the choice of the (dimensionless) temperature T* := a/lh 
and the (dimensionless) chemical potential p* := Pp. 
Within the present work we considered temperatures and 
chemical potentials in the ranges T* € [0.04,10] and 
p* G [—8,0], respectively. Box sizes were varied in the 
range L G {8 ,..., 19}. 

From an analysis of the equation of state, which are 
investigated as isotherms in q-p*-space, and of the struc¬ 
ture in terms of the pair distribution functions the phase 
diagram of the LRPM for various sets of values of the 
decay constant a and of the cut-off radius Tout are in¬ 
ferred. The following subsections explain the procedures 
to estimate the locations of the binodals of the first-order 
gas/liquid phase transition and of the so-called A-line 
of the continuous phase transition between the charge- 
ordered and the charge-disordered phase. 

C. Estimation of the binodals 

Estimates of the binodals of the first-order gas/liquid 
phase transition are obtained by determining the equa¬ 
tion of state t]{T*,p*). Figure 1 displays isotherms 
T* G {1.25,1.43,1.67} for L = 12,a = 0,rcut/cr = 1 us¬ 
ing the 1-norm (p = g = 1) in Eq. (2). For the lowest tem¬ 
perature T* = 1.25 a discontinuity in the packing frac¬ 
tion q{p*) indicates a first-order gas/liquid phase transi¬ 
tion, whereas for the highest temperature T* = 1.67 no 



r/o- r/cr 


FIG. 2. Pair distribution functions for L = 12, a = 
0, rcut/u = 1 using the 1-norm (p = g = 1) in Eq. (2) at 
temperature T* = 2 and packing fraction q = 0.37 (a) and 
q = 0.90 (b). The red curves represent the pair distribution 
function for equally charged ions, while the blue curves rep¬ 
resent those of oppositely charged ions. Panel (a) shows the 
structure of the charge-disordered phase, which is character¬ 
ized by a rapid vanishing of the spatial correlations. Panel 
(b) displays the case of the charge-ordered phase, where cor¬ 
relations are long-ranged due to a shell-wise alternating dis¬ 
tribution of cations and anions. 


gas/liquid phase transition occurs. An estimate of the 
critical temperature T* is obtained as a temperature for 
which q{p*) is discontinuous for T* < T* and continuous 
for T* > T*. This simple method is sufficiently precise 
for the present purpose. 

In order to estimate the binodal packing fractions 
qi{T*) and q 2 {T*) of gas/liquid coexistence at a tem¬ 
perature T* < T* the phase transition is first located 
within a range of the chemical potentials p* which is 
as narrow as possible and outside which the well-known 
hysteresis effect of first-order phase transitions does not 
occur. For T* = 1.25 in Fig. 1, e.g., the inset leads 
to the range p* G [—2.33,-2.28]. Our estimate of the 
chemical potential p*oe^{T*) at coexistence is the mid¬ 
point of the p*-interval of steepest increase of the pack¬ 
ing fraction q{p*). For the isotherm depicted in the in¬ 
set of Fig. 1, this choice yields /i*oex(^*) ~ —2.305 for 
T* = 1.25. Estimates of the binodal packing fractions 
771 (T*) and 772 (T*) are obtained by linearly extrapolat¬ 
ing the hysteresis-free parts of q{p*) to p* = pt^e^iT*) 
from below and from above, respectively. For Fig. 1, e.g., 
these estimates are qi{T*) « 0.27 and q 2 {T*) ~ 0.73 for 
T* = 1.25. Again, this simple method is sufficiently pre¬ 
cise, since, for the purpose of the present work, the loca¬ 
tion of the binodals is required to be known only semi- 
quantitatively but not to the highest technically feasible 
accuracy. More sophisticated and thus precise methods, 
such as umbrella sampling techniques or methods based 
on finite-size scaling and histogram reweighting [33, 34], 
are not required here. 
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FIG. 3. Order parameter |0| as a function of the chemical 
potential fj,* for decay constant a = 0 , cut-off radius rcut/o' = 
1 using the 1-norm (p = g = 1) in Eq. (2), and temperature 
T* = 1.67 for box sizes L € {12,18}. It signals the second- 
order phase transition between the charge-disordered (CD) 
phase, where 101 « 0, and the charge-ordered (CO) phase, 
where |0| 96 0. In order to account for finite-size effects, the 
convention is adopted that the phase transition occurs at | 0 | = 
0.03, which is located at p* « —1.65 here. 


D. Estimation of the A-line 

The structure of the charge-ordered and charge- 
disordered phase can be observed by analyzing the pair 
distribution functions gij{r) (see Eq. (4)). Figure 2 
shows the respective structures of both phases for the 
case r* = 2,L = 12,a = Ojrcut/o" = 1 using the 1-norm 
(p = g = 1) in Eq. (2). For the charge-disordered phase, 
spatial correlations vanish rapidly, while for the charge- 
ordered phase one observes long-ranged correlations due 
to a shell-wise alternating assembly of cations and anions. 

In order to locate the charge-ordered/charge- 
disordered phase transition, i.e., the A-line, in the phase 
diagram, the so-called staggered order parameter 



is considered. \(j>\ is positive in the charge-ordered phase 
and it vanishes in the charge-disordered phase. In or¬ 
der to account for finite-size effects, the convention is 
adopted that a thermodynamic state belongs to the 
charge-disordered phase iff |0| < 0.03. Obviously, this 
convention is largely arbitrary, but it allows a sufficiently 
precise localization of the A-line. Figure 3 displays \(j)\ 
as a function of the chemical potential p* for decay con¬ 
stant a = 0, cut-off radius Tcut/o = 1 and temperature 
T* = 1.67 for box sizes L G (12,18). For this example, 
the second-order phase transition, which is expected to 
belong to the Ising universality class [35], is located at 
fi* w -1.65. 



FIG. 4. Critical temperature T* as a function of the cut-off 
radius rent for sharp cut-off truncation schemes, where the 
decay constant a = 0 is used in Eq. (1). Data points labeled 
by T-norm’ and ’2-norm’ are obtained by truncation of the 
Coulomb potential at distances r^ut/o < {L — l)/2, where 
all distances are measured in the respective norms, i.e., p — q 
(see Sec. Ill A). Moreover, the values within the minimum im¬ 
age approximation (’MIA’) are displayed, which corresponds 
to p = 2, (jr = oo,rcut/o = {L — l)/2 (see Sec. IIIB). By 
comparison with the critical temperature 0^ G [0.14, 0.16] for 
the long-ranged Coulomb potential (rcut = 00 ) 9ne observes 
that the MIA results match the long-range value very well, 
whereas the 1 -norm and the 2 -norm results largely overesti¬ 
mate the critical temperature. 


III. RESULTS AND DISCUSSION 
A. Sharp cut-off schemes 

In this subsection a sharp cutting off of the Coulomb 
potential, i.e., a = 0 in Eq. (1), is considered with cut-off 
radii Tcut/o < {L — l)/2. First, the discussion is con¬ 
cerned with the dependence of the critical temperature 
T* on the cut-off radius rcut, which is obtained following 
the procedure described in Subsec. IIC. 

Figure 4 compares T* using the I-norm (p = q = 1) 
and the 2-norm {p = q = 2) in Eq. (2) with the crit¬ 
ical temperature 0* G [0.14,0.16] of the long-ranged 
Coulomb system, i.e., for rcut = 00 [30, 31, 36]. Distinct 
oscillations of T* for both 1-norm and 2-norm can be ob¬ 
served, where the amplitude of the latter is much smaller 
than that of the former. For both metrics, no signs of con¬ 
vergence of T* towards the value 0* of the long-ranged 
potential are observable within the considered range of 
cut-off radii rcut- Moreover, the amplitude of the oscilla¬ 
tions appears to even increase for the 1-norm. This odd- 
even-dependence likewise occurs within the structure for 
sufficiently large packing fractions: While for odd val¬ 
ues of rcut/o" a charge-ordered phase is realized at large 
packing fractions, no charge-ordering is observed for even 
values of rcut/c. 

In order to understand the presence or absence of a 
charge-ordered phase in the case of odd or even values 
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of ?'cut/o’, respectively, consider the particular case of the 
1-norm metric {p = q = 1 in Eq. (2)) and a perfectly 
charge-ordered configuration Cco- The energy of such a 
configuration is given by (see Eq. (2)) 

t 3 c_i y 

^HKco] = 2^ X! - i -^shell(*cr) 

/ 18 38 66 102 \ 

( 6 ) 

where 14heii(*o’) = -|- 2 is the number of sites in 

the i-th shell of 1-norm distance ia around a site. For 
?'cut/<7 = 1 only nearest neighbors interact and due to 
PH[(co] < 0 it is energetically favorable for the ions to 
arrange alternatingly for sufficiently large rj. However, 
for rent/o' = 2 the colons in the second shell overcompen¬ 
sate the contribution of the counterions in the first shell, 
such that PH[(^co] > 0, which renders the charge-ordered 
phase unfavorable as compared to the charge-disordered 
phase. This pattern repeats upon increasing Tout with 
PH[Cgo] < 0 for odd and /3H[(co] > 0 for even values of 
rcut/o. Since the size of the i-th shell grows oc and the 
strength of the interaction with the central ion decreases 
oc 1/i, the total energy PH[(co] oscillates with an am¬ 
plitude which increases oc i. This growing stability of 
the charge-ordered phase for odd values of Veut/o causes 
the increase in the critical temperature T* for odd-valued 
cut-off radii. Furthermore, for all cases studied here, it is 
found that the total energies within the charge-ordered 
phases for odd-valued cut-off distances are much lower 
than for the long-ranged Coulomb potential and therefore 
that the stability of the latter’s charge-ordered phase is 
much weaker. This leads to a significantly lower critical 
temperature compared to those, obtained using the sharp 
cut-off scheme. The absence of a charge-ordered phase 
in the case of even-valued cut-off radii and the increas¬ 
ing deviation of the critical temperature T* from that of 
the conventional LRPM with long-ranged Coulomb po¬ 
tential, 0* S [0.14, 0.16] in the case of odd-valued cut-off 
radii leads to the conclusion of some major defect of the 
sharp cut-off scheme. 

Calculations using the 2-norm lead to qualitatively 
similar results. The quantitative differences for Ccut /o > 
1 to the case of a 1-norm metric can be attributed to the 
occurrence of non-integer-valued distances and to the fact 
that the shells of constant distance from a central site 
are spheres and not octahedrons. Note, that for the case 
J'cut /o = 1 the distinction between 1-norm and 2-norm is 
irrelevant. Remarkably, the critical temperatures T* for 
odd values of Veut/o are considerably smaller than those 
for the 1-norm metric, whereas for even values of rcut/o 
the differences are very small. 

Figure 5 displays the phase diagrams for cut-off radii 
’’’cut/o € {1,3,5} (1-norm, p = q = l), which have been 
determined using the methods presented in Secs. IIC and 
IID. Close to the critical point, the binodals follow a 



FIG. 5. Phase diagrams for 1-norm metric (p = q = 1 in 
Eq. (2)), decay constant a = 0, box size L = 12, and cut-off 
radiircut/cr € {1,3,5} compared with that of the conventional 
long-ranged LRPM (ccut = oo [31]). All these phase diagrams 
exhibit the same topology of a charge-ordered and a charge- 
disordered phase separated by a first-order phase transition 
below and a continuous phase transition above a tricritical 
point (•). However, for the cases with finite values of the cut¬ 
off radius rcut, the tricritical temperature is about an order 
of magnitude too high and it increases upon increasing Ccut- 


straight line, which corresponds to the critical exponent 
/3 = 1. Due to this observation and since the A-line termi¬ 
nates at the critical point, this critical point actually is a 
tricritical point (tagged by a black dot). The long-ranged 
LRPM qualitatively shows the same phase diagram, in¬ 
cluding tricriticality. In accordance with the results of 
Fig. 4, the tricritical point moves upwards to higher tem¬ 
peratures for larger odd values of rcut/o', however, the 
topology of the phase diagram is not affected. It has al¬ 
ready been mentioned that for even values of the cut-off 
radius rcut/o no charge-ordered phase is present so that 
a phase separation occurs between a charge-disordered 
gas and a charge-disordered liquid, and the phase dia¬ 
gram exhibits an ordinary critical point and no A-line. 
Since the phase diagrams for even-valued cut-off radii 
are qualitatively different from that of the long-ranged 
LRPM they are not shown in Fig. 5. 


B. Minimum image approximation (MIA) 

When discussing the sharp cut-off schemes, i.e., with 
decay constant a = 0 in Eq. (1), in the previous Sub¬ 
sec. Ill A, cut-off radii rcut/o < {L — I)/2 have been 
assumed throughout. The phase diagram for odd val¬ 
ues of the cut-off radius rcut/o < {L — 1)/2 turned out 
to be qualitatively identical to the phase diagram of the 
conventional long-ranged LRPM, but the tricritical point 
is located at a too high temperature T*, which even in¬ 
creases upon increasing rcut- Here, another sharp cut-off 
scheme is discussed which, in the notation of Eq. (I), can 
be specified by a = O,^ = oo,rcut/o = [L — I)/2. This 
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FIG. 6. Comparison of the equations of state and of the pair 
distribution functions for odd and even box sizes L. The green 
and blue dots in the main plot correspond to results within 
MIA {p — 2, q — oo, Q = 0, rcut/cr = {L — l)/2 in Eq. (1)) and 
the black dots are obtained by means of the Ewald method. 
The same odd-even-effect of L in the high-density regime is 
encountered within both methods. Furthermore, the high- 
density phase displays decorrelations in the case of odd values 
of L, whereas the charge ordering is long-ranged in the case of 
even values of L. The occurrence of decorrelations in the high- 
density phase for odd-valued L can be understood in terms of 
equally charged ions at the rim of the simulation box, which 
interact with each other due to periodic boundary conditions. 


cut-off scheme is equivalent to the well-known minimum 
image approximation (MIA) [37, 38], which disregards 
all contributions to the interaction energy, which do not 
correspond to the minimum distance between two ions or 
their periodic images. For a cubic simulation box, it can 
be interpreted such that only those interactions within 
a cut-off distance ||r||oo < »"cut are taken into account, 
where the cut-off radius rcut has to be chosen such that 
the interaction range equals a cube of volume with 
L = 2rcutlcy + 1- Note, that the distances, which deter¬ 
mine interaction potential, are measured in the p-norm 
(see Eq. (1)), which is chosen as p G {1,2}. 

Figure 4 also displays the critical temperatures T* 
within MIA for p = 2, which are in excellent agree¬ 
ment with the critical temperature 0* € [0.14,0.16] of 
the long-ranged LRPM. Moreover, the MIA results in 
Fig. 4 can be subdivided into integer and half-integer 
values of rcut/c, which corresponds respectively to odd 
and even values of the box size L = 2rcutler -|-1. A close 
look at Fig. 4 reveals that the values of T* within MIA 
for odd values of L are slightly below those for even val¬ 
ues of L, which indicates a certain odd-even effect of L 
on T*. Whereas the odd-even effect of L on T* is rather 
weak. Fig. 6 shows for the examples L = 13 and L = 14 
that there is a significant odd-even effect of L on the 
equation of state and on the structure at large packing 
fractions. The reason for mentioning the odd-even effects 
of L within MIA here is that exactly the same odd-even 
effects of L occur for the long-ranged LRPM using the 
full Ewald method, which are shown in Fig. 6, too. In 


other words, MIA is such a good approximation of the full 
Ewald method that it exhibits even the same odd-even 
effects. 

The insets in Fig. 6, which show the pair distribution 
functions at the indicated thermodynamic states, reveal 
a distinct charge-ordering behavior in the high-density 
regime for even box sizes L, whereas for odd L a decaying 
oscillatory behavior of the pair distribution functions can 
be observed. Although the latter behavior does not cor¬ 
respond to a genuine charge-ordered phase, it is clearly 
different from the structure in the dilute regime, where 
correlations decay rapidly. What is the reason for this 
qualitatively different high-density structure for even and 
odd box sizes L? In the limit of perfect charge-ordering at 
sufficiently high densities, the periodic boundary condi¬ 
tions lead to destabilizing contributions to the interaction 
energy for odd values of L, since in this case there are 
equally charged ions located at the rim of the simulation 
box which, due to the periodicity, are nearest neighbors. 
This effect does not occur for even box sizes L. Since 
this is obviously an effect related to the ‘surface’ of the 
simulation box, it decreases with increasing box size L. 
The same odd-even dependence of the equation of state 
and of the structure is observable within the long-ranged 
LRPM. As a consequence of the above argumentation, 
the results for even box sizes L, whether using MIA or 
the full Ewald method, can be expected to be more ac¬ 
curate, than for odd box sizes L. 

It is striking that isotherms within MIA sample the 
long-range limit rcut oo very well for dilute as well as 
for dense systems. This indicates that the Hamiltonian 
of an ionic fluid with long-ranged Coulomb interaction 
can be reliably approximated by Eq. (2) with p = 2,q = 
oo,a = 0, rcut/CT = {L — l)/2 in Eq. (1). In fact, for 
L = 14, there is almost no difference between MIA and 
the full Ewald method with respect to the equations of 
state (see Fig. 6). 

This quantitative agreement can be understood as fol¬ 
lows: The contribution to the energy due to the long- 
ranged Coulomb potential beyond MIA involves the elec¬ 
trostatic interaction of the simulation box with its im¬ 
ages, which, due to the charge neutrality of the simula¬ 
tion box, decays at least as a dipole-dipole interaction, 
i.e., cx 1/r^. However, numerical calculations reveal that 
the dipole-dipole contribution vanishes and that the de¬ 
cay is actually oc 1 /r^, which renders the long range con¬ 
tributions beyond MIA absolutely convergent and, for 
sufficiently large L, small. The same argument does not 
apply to the sharp cut-off potentials of Subsec. HI A, 
since no charge neutrality inside spheres of radius rcut 
within the underlying norm is guaranteed so that the long 
range contributions beyond those due to the sharp cut¬ 
off potential are due to an effective monopole-monopole 
interaction, i.e., oc 1/r, and, hence, typically not small. 

Figure 7 displays the phase diagram within MIA for 
feut/er = 4.5; those for rentier G {6,6.5,9} obtained sim¬ 
ilarly are not shown here. For all considered cases of 
Tout, the phase diagrams are qualitatively and quantita- 
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FIG. 7. Comparison of the phase diagrams of the LRPM 
within the sharp cut-off scheme p = q = l,aa = Ojrcut/o" = 
5 < (L — l)/2 in Eq. (1), the MIA {p = 2, q — oo,aa — 
0, rcut/cr = 4.5 = (L — l)/2), the smooth cut-off scheme 
{p = q = 2,aa = 0.8, rcut/cr = 5 < {L — l)/2), and the 
long-ranged Coulombic interaction {aa = 0, rent = c)o). The 
agreement of the MIA and the smooth cut-off scheme with the 
case of the long-ranged interaction is excellent, whereas there 
are large quantitative deviations of the sharp cut-off scheme 
with rcut/o- <{L — l)/2. 


tively in accordance with those obtained by means of the 
corresponding full Ewald method. The important point 
is that the MIA can be interpreted as a sharp cut-off 
truncation scheme of the Coulomb potential. The results 
discussed above indicate that the long-range character 
of the electrostatic interaction may not be necessary for 
the thermal and structural properties of ionic fluids, be¬ 
cause there is the possibility that the same properties can 
be generated by suitable short-ranged interactions. A 
comparison of the MIA, where rcut!o' = [L — l)/2, with 
the sharp cut-off schemes of the previous Subsec. Ill A, 
where Veut/o < {L — l)/2, reveals that the ability of a 
sharp cut-off scheme to mimic the properties of the long- 
ranged LRPM depends delicately on the relation between 
the cut-off radius rcut and the size of the simulation box 
L. However, the simulation box L is not a physical but 
rather a technical parameter, and the properties of the 
LRPM within reasonable sharp cut-off schemes Eq. (1) 
with a = 0 should be independent of (large values of) 
L. Consequently, there are no sharp cut-off schemes 
which quantitatively reproduce the thermal and struc¬ 
tural properties of the long-ranged LRPM and which, at 
the same time, correspond to physically acceptable un¬ 
derlying short-ranged interactions. 


C. Smooth cut-off potentials 

It turned out in the previous Subsec. IIIB that 
there are sharp cut-off schemes, i.e., with a = 0 in 
Eq. (1), which, although based on a short-ranged po¬ 
tential Eq. (1), lead to quantitatively identical thermal 



aa aa 


FIG. 8. Critical temperature T* and critical packing fraction 
Pc of the LRPM with smooth cut-off potentials Eq. (1) with 
p — q — 2, L = 12, rcut Io G {1)5} as functions of the decay 
constant aa. For Veutja = 5 the critical temperature T* ex¬ 
hibits a plateau around aa « 0.8 which is close to the critical 
temperature 0* of the long-ranged Coulombic system. As in 
Figs. 4, 5, and 7, the critical temperature T* within sharp 
cut-off schemes [aa = 0) is an order of magnitude to high. 
The critical packing fraction pc ~ 0.5 is rather independent 
of the decay constant a and the cut-off radius rcut. 


and structural properties as the conventional LRPM with 
the long-raged Coulomb potential. However, these sharp 
cut-off schemes do not correspond to physical systems, 
because one has to tune the cut-off radius rcut accord¬ 
ing to the unphysical simulation box size L, whereas for 
physically acceptable short-ranged potentials rcut should 
be independent of L. Here, potentials Eq. (1) with a > 0 
are considered, which correspond to a smooth cutting 
off of the Coulomb potential on the length scale l/ct. 
Note, that this type of functions occurs within the Ewald 
method as a result of splitting the total electrostatic po¬ 
tential into a short-ranged contribution, which leads to a 
sum in real space, and a long-ranged contribution, which 
leads to a sum in reciprocal space. Within the Ewald 
method, a is adjusted such that, on the one hand, the 
error due to the unavoidable truncation of the sums in 
real and reciprocal space are sufficiently small and that, 
on the other hand, the computational effort is accept¬ 
able. However, here a controls the decay length of a 
genuine short-ranged interaction, which can take, in prin¬ 
ciple, any (positive) value. Eor a —>■ 0 the sharp cut-off 
truncation schemes discussed in Sec. HI A are obtained. 

Figure 8 displays the critical temperature T* and the 
critical packing fraction pc of the LRPM with smooth 
cut-off potentials Eq. (1) with p = q = 2, r cut jo S {1,5} 
as functions of the decay constant ao. For rcut jo = 5 a 
plateau of T* is present in the range ao G [0.4,1.0], where 
the critical temperature T* is quantitatively equivalent 
to the critical temperature 0* of the long-ranged LRPM 
(see Fig. 8(a)). For larger values of ao the decay length 
of the interaction Eq. (1) becomes too short, so that even 
nearest neighbors barely interact with each other. Hence, 
for ao oo, one obtains the ideal gas limit and conse¬ 
quently T* —>^0. For rcut/o = 1 all contributions beyond 
the nearest neighbors are neglected, and no plateau is 









r/cr 


FIG. 9. Accumulated charge Q{r) (see Eq. (7)) up to a dis¬ 
tance r in 1-norm around a positively charged central ion. 
The Stillinger-Lovett sum rule [13-16] requires Q{r) to van¬ 
ish in the limit of large distances r —>■ oo. This condition ap¬ 
pears to be fulfilled within the charge-disordered phase (here 
at packing fraction 77 « 0.3) for the case of the smooth cut¬ 
off potential aa = 0.7 (p = q = 2,rcut/o' = 5, L = 12) as 
well as for the MIA {aa = 0,p = 2,q — oo,rcut/a = 6.5 = 
{L — 1)I2,L = 14), but not for the sharp cut-off potential 
acr = 0 (p = g = 2, r^utla = 1, L = 12). 

observable at all. However, the critical packing fraction 
rjc ~ 0.5 appears to be independent of the decay constant 
a and the cut-off radius Tcut- In addition to the coinci¬ 
dence of the critical points for the appropriate choice of 
aa, a charge-ordered phase at large packing fractions and 
a charge-disordered phase at small packing fractions can 
be found. Ultimately, the phase diagrams of the smooth 
cut-off scheme are in quantitative agreement with the 
long-ranged limit, too (see Fig. 7). 

The long-range character of the Coulomb interaction 
is well-known to generate the so-called Stillinger-Lovett 
perfect screening property of ionic fluids [13-18]. In the 
present context of an LRPM it corresponds to the accu¬ 
mulated charge 

Q(r):=l-|-^ ^ [5e,e(r0 “ 5e.e(i’0] (7) 

r'ev 

l<l|r'lli<’' 

up to a 1 -norm distance r around a positively charged 
central ion to vanish Q{r) —>■ 0 in the limit r ^ 00 . 
Figure 9 compares the accumulated charge Q{r) for the 
sharp cut-off potential with aa = 0-,rcut/(^ = liT = 9 = 
2, the MIA with aa = 0, rcut/o" = 6.5, p = 2,q = 00 and 
the smooth cut-off potential aa = 0.7, rcut/CT = 5,p = 
q = 2 within the charge-disordered phase at packing 
fraction 77 ss 0.3. For the sharp cut-off potential, Q{r) 
appears to converge towards a non-vanishing (positive) 
value for large radii r, which corresponds to an imperfect 
screening of the positive central charge. However, perfect 
screening occurs within the MIA and the smooth cut-off 
potential. This is an interesting finding since, in contrast 
to the long-range Coulomb potential, the perfect screen¬ 


ing property is not necessarily fulfilled in systems with 
short-ranged interactions. 

It can be shown for smoothly cut-off potentials Eq. (1) 
with Teuton 2 > 1 that the charge-charge pair correlation 
function /i^^(r) := 2 (p®_®(r) - p 0 ,e(r)), which is re¬ 
lated to the accumulated charge Q{r) in Eq. (7), decays 
asymptotically on the length scale of the Debye length 
1/k with (kc)^ = Anp/T* iff k ^ 2a. Since a > 0 is a 
constant, which is chosen to match the phase diagram of 
the long-ranged Coulomb interaction (see Figs. 7 and 8 ) 
and which does not depend on the thermodynamic state 
( 77 , T*), the condition k 3 > 2a > 0 will not be fulfilled at 
very low packing fractions 77 and/or high temperatures 
T*, i.e., the Debye-Hiickel limit is not recovered within 
smooth cut-off schemes. However, we do not consider this 
a major defect for two reasons: On the one hand, the re¬ 
maining parts of the phase diagram outside the region of 
the Debye-Hiickel limit, particularly in the range of high 
densities, e.g., close to the critical point {rjc ss 0.5), are 
reproduced quantitatively (see Figs. 7 and 8 ). On the 
other hand, the crossover, where the condition 7 t 3 > 2 a 
begins to be violated, can be shifted to arbitrarily small 
values of k by decreasing a, which can be achieved by 
accordingly increasing rcut (see Fig. 8 (a)). 

Concerning the phase behavior and the structure, the 
LRPM with the smooth cut-off potentials considered 
here is qualitatively equivalent to that with long-ranged 
Coulomb interactions. Moreover, by choosing an appro¬ 
priate decay constant a, the short-ranged smooth cut-off 
potential LRPM becomes even quantitatively equivalent 
to the conventional long-ranged LRPM. The essential dif¬ 
ference to the sharp cut-off schemes considered in Sub¬ 
secs. HI A and HIB is that these statements are indepen¬ 
dent of the choice of the actual (large) value of the simu¬ 
lation box L, i.e., the smooth cut-off potentials discussed 
here are candidates of physically meaningful short-ranged 
interactions. Consequently, the long-range character of 
the electrostatic potentials is of minor importance for the 
thermal and structural properties of ionic fluids, since, 
as has been demonstrated here, short-ranged interaction 
potentials do exist, which lead to the same thermal and 
structural properties. 


IV. CONCLUSIONS AND SUMMARY 

In order to elucidate the role of the long-range char¬ 
acter of the Coulomb potential for the thermodynamic 
and the structural properties of ionic fluids, the lattice 
restricted primitive model (LRPM) with various trun¬ 
cated Coulomb-like interaction potentials Eq. (1) has 
been studied in this work by means of grandcanonical 
Monte Carlo simulations. 

The simplest approach of sharply truncating the 
Coulomb potential at a certain cut-off radius rcut 
(Sec. HI A) turned out to be not appropriate to reproduce 
the properties of the long-ranged LRPM, either qual¬ 
itatively or quantitatively. When measuring distances 
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within the 1-norm, e.g., no charge-ordered phase occurs 
for even-valued cut-off radii rcut/o", i.e., the topology of 
phase diagram differs from that of the conventional long- 
ranged LRPM. On the other hand, odd-valued cut-off 
radii lead to qualitatively the same phase diagrams as 
the conventional long-ranged LRPM (Fig. 5), but there 
is a large quantitative difference, e.g., of the value of the 
critical temperature T* (Fig. 4). 

In contrast, the well-known minimum image approxi¬ 
mation (MIA), which may be viewed as a sharp cut-off 
scheme with the cut-off radius Tout being related to the 
side length of the simulation box L, leads to quantitative 
agreement of the phase behavior and the structure of the 
LRPM with that of the long-ranged LRPM (Fig. 7). 

However, due to the relation of the cut-off radius with 
the unphysical parameter L, the MIA does not corre¬ 
spond to a physically acceptable, i.e., L-independent. 


short-ranged potential. Quantitative agreement of the 
phase behavior and the structure of short-ranged and 
long-ranged LRPMs is achieved by smoothly cutting off 
the Coulomb interaction, which is expressed by a de¬ 
cay constant a, corresponding to a decay length of 1/a 
(Fig. 8). Interestingly, the LRPM with the smoothly 
cut off Coulomb potential possesses, although it is based 
on short-ranged potentials, the Stillinger-Lovett perfect 
screening property (Fig. 9). 

By showing that the characteristic phase behavior and 
structure of ionic fluids can also be found in systems with 
short-ranged potentials, one can conclude that the deci¬ 
sive property of the electrostatic potential in ionic fluids 
is not the long range but rather the valency dependence. 
As this conclusion is drawn from a study of the LRPM, 
the natural question arises whether it holds for off-lattice 
models, too. This task is left for future investigations. 
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